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tG ' A simple variational Lagrangian is proposed for the time develop- 

ment of an arbitrary density matrix, employing the " factorization" of 
the density. Only the "kinetic energy" appears in the Lagrangian. 
The formalism applies to pure and mixed state cases, the Navier- 
1^ ' Stokes equations of hydrodynamics, transport theory, etc. It recap- 

tures the Least Dissipation Function condition of Rayleigh-Onsager 
and in practical applications is flexible. The variational proposal 

\^^ , is tested on a two level system interacting that is subject, in one in- 

^^ ' stance, to an interaction with a single oscillator and, in another, that 

\l , evolves in a dissipative mode. 

o 

^ '■ 1 Introduction 

*f^- Several basic aspects of stochastic dynamics remain controversial. (A critical 

Q_i, "state of art" update in [1] shows this). This situation contrasts with most 

physical theories, where the problems that arise are in the application of 
consensually accepted principles. It can perhaps be argued, that the lack of 
an agreed-upon variational formulation of stochastic processes is at the root 
d I of the problem. As a remedy of this situation, this article suggests a new 

variational functional which is to be minimized and whose minimum is the 
true density matrix. 

To be sure, in the past several principles of extrema have been proposed; 
these include Gauss's Least Constraint Principle [2, 3], the "least dissipation 



function" [4, 5, 6], minimum entropy production-rate for steady states [7, 8] 
(see also [9] for its violation), minimal energy generation rate [10], minimal 
scattering integral ([11] - [13]), least velocity error during pathway [3], the 
Yasue Action for stochastic mechanics [14], a formulation involving a poten- 
tial [15], and again, recently, maximum entropy production [16]. To these 
may be added several variational methods applicable to classical (i.e., not 
quantal) systems, such as those appropriate for general non-linear problems 
[17], the "governing principle for Dissipative Processes" [18, 19] and a gen- 
eralized Hamiltonian principle [20] . Reviews of these and of other methods 
can be found in [21] - [24]. 

The present proposal for a variational procedure is based on the follow- 
ing new elements: (a) the factorization (to be discussed later in this paper) 
of the density matrix as introduced by Reznik [25] and utilized recently by 
Gheorghiu-Svirschevski [16], (b) a conventional Lagrangian similar to that 
used in Mechanics to obtain the motion of a point particle subject to an ex- 
ternal force, but in which the scalar potential is either absent or ignorable, 
(c) a vector potential that can be singular, without this having disastrous 
observational consequences and (d) an appropriate use of minimization pro- 
cedure, with origins going back to Gibbs [2]. The method covers a broad 
range of behaviors ("deterministic" and stochastic, quantal and classical, 
electronic transport, discrete and continuous, markovian and otherwise) and 
places in a new perspective certain aspects in currently employed theories of 
stochastic dynamics. Apart from these favorable (and a priori unexpected) 
features, a number of problems remain, which will have to be resolved by 
future efforts. 

A pure state sub-case, (as opposed to a density matrix that can de- 
scribe a statistical mixture of states) is the subject of a numerically worked 
out example in section 4.1 and is equivalent to the (linear) time dependent 
Schrodinger equation. For this several variational formulations are known in 
the literature ([26] -[31]). The inter-relation between these was investigated 
in [30], where they were shown to be frequently equivalent. For the pure 
state case our density matrix variational method reduces to the McLachlan 
formalism [27], in which the variation of the function is carried out with 
respect to the time-derivative only (while the function is kept fixed). More- 
over, we give a justification of this procedure for stochastic processes. 

Another application of the variational formalism, in section 4.2, includes 
a non-linear, dissipative (non-Hermitian) mechanism and exemplifies quan- 
tum jumps. 



2 The Variation 

The "factorized" form of the time dependent density matrix p ([25], [16]) 
reads in terms of the column and row vectors 7 and (its hermitian conjugate) 

7+ 

P = 77^ (1) 

The above condensed notation is not trivially simple, so we give in Appendix 
A a "tutorial" on the notation. We now propose an action S{T) (expressed 
in arbitrary units) this being the integral over time t with an arbitrary time 
end-point T of a Lagrangian C. It is the form of the Lagrangian that we 
seek: we propose that it has the "quadratic" form, as follows. 

S{T) = f dt C{t) 
Jo 

= f dtTr (i7+ - A+) • {-ij - A) (2) 

The variational equations based on the above action (these are the equations 
of motions for 7 and 7''') are given below in equation (7) . Tr is the trace 
over all "components" of 7 explained in Appendix A. Dots above symbols 
represent differentiation with time. 

The variational equations are obtained in the usual way by varying the 
action S{T) with respect to all components of the two vectors 7 and 7'''. 
Thus 



6,SiT) = j^ m-,m-if - a^^l + ^h(T) (3) 

and a similar expression for 5^+S{T). The last term, outside the integral, 
is the boundary term. It is assumed that the vectors 7 are fixed initially 
at t = 0, but not at any time later. (This will be discussed shortly.) If the 
boundary term can be made to vanish, then so will be the whole variation, 
since in the absence of a scalar potential, the variations ^^ ' and . \' 
vanish. This follows, since for the postulated form of the Lagrangian these 
variations contain as a factor one or the other of the expressions g. ' and 

..\' and these vanish due to the boundary variation. (In fact, —^ and 

-T-4p- also vanish due to the presence of the same factors, but these variations 
do not form a sufficiently general basis for the variation procedure, since the 
vector potential A may not be a function of all components of the 7-vector.) 
At this point the role of the boundary term is well worth reflecting 
upon. It is not present in, e.g., deterministic Mechanics, where the values 



of the variables are fixed at both end points. However, it is well known 
that the boundary term arises when the value of the variant quantity, i.e. 
7, is undetermined at a boundary. This is (physically) the case when a 
random force operates on the system. Thus, we are not allowed to neglect 
the boundary term. It is now a further "bonus" in the formalism, that 
the vanishing of the boundary-term variation does not interfere with the 
vanishing of the body terms- variation, (i.e., it is neither contradictory to 
it nor incomplete with respect to it), but is by virtue of our choice of the 
Lagrangian precisely identical with it. 

In summary, we have the variational equations 

^ = (4) 

and their complex conjugates 

5C{t) 

6j+ 



(5) 



and these make the action extremal (and an absolute minimum) also when 
the "vector potential""^ A is a function of 7 [32]. [To see this, note the remark 
after equation (7) ] This has the immediate consequence that in the expan- 
sion of the Lagrangian, shown in equation (2) , the following expression 
needs to be minimized: 

'DF' = 7+7 - i[7+A - 7A+] (6) 

This follows, since the A • A^ term is independent of time-derivatives and 
is not varied. We further note that (by the form of the vector potential) 
i times the square bracket is a real quantity. The quantity in the above 
equation is essentially of the form of Onsager and Machlup's Dissipation 
Function [the negative of eq. (4-25) in [5]. One recalls that their Dissipation 
Function is also minimized only with respect to the time derivatives of the 
variables just as in the procedure proposed in [27] and in the present work]. 
To bring our 'DF' precisely to the form of the Dissipation Function, we 
need to go from our variables 7 and 7^ by a constant linear transformation 



^We have named our frequently used quantity A the "vector potential", in analogy 
with the quantity that enters as a cross-term with the time derivative ("the particle 
velocity") in the Lagrangean of classical mechanics [32], or inside the square with the 
canonical momentum and in distinction from the scalar potential <j}, which appears in 
(non-relativistic) Lagrangeans as a separate term. 



(not necessarily a unitary one) to the variables a of [5]. One will then get, 
instead of the diagonal form 7^7, a non-diagonal form which defines the 
"generalized resistance matrix" Rij of [3, 6]. [Since the use of extensive 
variables is to be preferred to intensive ones (and 7 and 7"*" are of the latter 
type) the transformation should include the system size. Alternatively, the 
action integral may be premultiplied by a size-dependent scale factor.] 
Thus, the result of the variation are a set of equations: 

7 = iA 
7+ = -iA+ (7) 

When A and A"*" are functions of 7 and 7+, their (non-zero) derivatives 
come in with either 7 — iA or 7"'" -|- iA+ as factors and these factors vanish 
due to the above equations. 

[We can illustrate this in the case of two components, for which we write 
the Lagrangian in equation (1) as 

m = fih^fth (8) 

and /i = —^71 — Ai(7fc,7^) ^tc. with the 7 dependence explicitly put in 
the A's. Recall that equation (7) means that for the variational solution 
fi=fi=f2 = f2—^^^ ^11 times and that with this solution the minimal 
action is zero. 

Upon varying with respect to e.g., 71, the resulting Lagrange-Euler equa- 
tions is now in full detail 

-{d^,At)h + ^A+ - ftid^i^i) - id^iAt)f2 - ftid^^M) = (9) 

Since all the terms contain one of the /-factors or their time derivatives 
(which are necessarily also zero) the above equation is satisfied for the pro- 
posed variational solution. There may be other solutions, though! If for 
these not all /-s are zero (and therefore differ from the proposed solution), 
than the action (which consists of positive terms) is positive and larger than 
that for the solution given in equation (7) .] 

These are the equations of motion of the (independent) vector variables 
and can be regarded as having the status of the Langevin equations, or 
Hamilton's equations for the set of conjugate variables 7a and 7^^. The 
processes considered in this section comprise a purely Hamiltonian process 
Ah (namely, energy preserving, "elastic") as well as some other, dissipa- 
tive mechanisms. Thus, a Markovian scattering process represented by the 



symbol M^a (designating half the probability per unit time of a scattering 
event taking the system from state a to b) can be written as {Aout,Ain) to 
separate scattering out of and into a given state. We also add a stochastic, 
random process arising from, e.g., an external source, as Ar. Some other 
type of processes will be considered below in section 5. 
For the first two processes we have 

Ah = -H-f 

A+ = -7+F+ = -7+i7 

iAout)a = i^MbaJa 
b 

b 

b 

{Al)a = iha)-H^J2Mab7blt) (10) 

b 

When we add to these the random force, we obtain in addition 

iAr)a = -ifala (^+)a = llafa (H) 

fa represent the components of the random time-dependent force with zero 
mean and a finite self-correlation. (iV is the number of states in the ensemble, 
see appendix A) 

The vector potential Ain is singular. However, singularities in vector 
potentials are well known (as, e.g., in those for solenoidal or monopole fields). 
To cancel these singularities we shall follow the procedure of Reznik [25] 
and Gheorghiu-Svirshevski [16], who multiply 7"^ into 7, 7 into 7^^, add and 
obtain the (master) equations for the density matrix. 

After substituting the quantities from equations (7,10,11) we obtain by 
this procedure for a diagonal element (say) aa of the density matrix 

Paa = (77^ + 77^)aa = i{[p, H])aa " '^Paa ^ M^a + 2 ^ MabPbb + '^Paafa 

b b 

(12) 

When one writes out the equation, similar to equation (12) , for the 

time derivatives of the off-diagonal matrix elements Pab{o- 7^ b), one finds 

singularities in them, due to the above mentioned singularities in the vector 



potentials. For a macroscopic system these singularities cancel, when one 
takes into account the phase decoherence between different states of macro- 
scopic bodies. (The subject of microscopic to macroscopic transitions does 
not belong here. It was studied by various methods and has summaries in 
e.g. [33, 34].) 

3 Potential Fluid Dynamics 

An interesting application of the preceding complex factor-density formal- 
ism is for the well known potential flow (namely, fluid dynamics without 
vorticity) as presented in many hydrodynamic text books, e.g., [35]. If a 
flow satisfies the condition of zero vorticity, i.e. the velocity field v is such 
that V X v = 0, then there exists a function (p satisfying v = V0. 

In this section we answer the question: what form of the vector potential 
A appearing in equation (2) will ensure, that upon variation of the action 
containing these vector potentials we shall obtain precisely the well known 
equations of potential flow hydrodynamics. These equations are: 

^ + V.(pW) = (13) 

f = -^{V^f-h-<^-.V'cP (14) 

In these equations the physical meaning of the quantities is that p is the 
mass density, h is the specific enthalpy, u is the viscosity coefficient and 
<I> is some function representing the potential of an external force acting 
on the fiuid. The first of these equations is the continuity equation, while 
the second is a modified Bernoulli's equation which takes into account some 
viscous effects (A full viscous fiow is of course not a potential fiow and 
contains vorticity) . p and (p play the roles of the squared-amplitude and of 
the phase angle, respectively. Both are real quantities. 

The final results for the desired vector potentials and their complex con- 
jugates are shown, below, in equation (23) and equation (24) . To obtain 
them, we first express the variational variables 7 and 7+ that we have used 
so far in terms of the physical variables (p, (p) . The variation will now be 
carried out with respect to the latter variables. The transformation is 

7 = ^pe'^ (15) 

7+ = ^pe-''f> (16) 



Though ah variables are now functions of the positions, and are thus con- 
tinuous variables, we shall label them, as before, by the subscripts a, etc. 
The following relations (with no summations over repeated symbols) arise 
simply from the inverse transformation: 

0a = {ial'"a-lait)l'^ipaa (17) 

V0a = (V7a7+a - 7aV7+ J/2ip„, (18) 

-{{^lal^a? - (7aV7+j2)/2^pL (19) 

Vpaa = (V7a7+a + 7aV7+ J (20) 

Using these, we first rewrite equation (13) as 

Paa = {ial^a + lalt) = -paa^'^'Pa " "^ Paa ' "^ 4>a 

^ Railan-^a) (21) 

Ra being a well-defined, real function of the variational variables and of their 
first and second spatial derivatives (but independent of the time-derivatives) . 
Likewise, one obtains as the rate equation for the phase, equation (14) , as 

^a = -{ial^a-lait)/'^ipaa 

= -^{V<pf -h-<^-uV^4> 

^ Naila.l^a) (22) 

in which Na has properties similar to Ra- 

We can solve for the two quantities, 7a7'''a a^id '^ait from the preceding 
two equations, and then divide by 7"'"q and 7a, respectively, to obtain the 
time derivatives. However, by equation (2) the time derivatives are just the 
vector potentials. Thus we finally obtain 

7a = {Ra + 2ipaaNa)^r^ = iAa (23) 

27+a 

and the complex conjugates 

7a+ = ^{Ra- 2ipaaNa) = -iA+ (24) 

27a 

We have thus found the vector potentials which have to be inserted in the 
action, so as to yield variationally the hydrodynamic equations, equation 



(13) and equation (14) . It is evident that the complex representation is a 
natural way to obtain variationally equations of motion for two such dissim- 
ilar quantities as amplitude and velocity. The physical extrema are certainly 
global minima (although the functional may have additional minima). De- 
tailed applications will be undertaken in the future including the problem 
of a general viscous flow. 

4 Applications 

4.1 A Periodically Varying Hamiltonian 

We shall now apply the proposed variational procedure to yield, in one case 
exactly and in another case approximately, the solution for a (hermitian) 
Hamiltonian that has a periodic variation. The cases chosen are such that 
analytical solutions are known exactly ([37]- [39]), so that we can compare to 
them the variational solutions to be obtained here. Specifically, we consider 
the time development of a doublet subject to a Schrodinger equation whose 
Hamiltonian in a doublet representation is 

\ stn[ujt) cos[ujt) I 

Here to is the angular frequency of an external disturbance. The eigenvalues 
of equation (25) are y and — y. If G > 0, then in the ground state the 
amplitude in the upper component ( j^ ) in equation (25) is 



C^ = cos{Kt) cos{ujt/2) + {uj/2K) sm{Kt) sm{ujt/2) 

+ i{G/2K)sm{Kt)cos{ujt/2) (26) 

with 

K = OWG^+uj^ (27) 

The amplitude C„ of lower component ( '^' ) in the ground state has a similar 
form, which we shall not bother to write out. For the variational procedure 
we postulate a superposition of complex circular functions 

C7«= Y, ^™exp^(™-5M (28) 

m=-M,M+l 

and similarly for CL The complex coefficients A^^ and ^J„ are determined 
by minimization of the variational action, equation (2) , subject to the nor- 
malization condition that the sum of the absolute squares of the coefficients 







Figure 1: Time-dependence of the real part of the amphtude in the upper 
electronic component state. Parameter values G = 3, a; = 1. Full line: 
variationally obtained state. Broken line: the exact, algebraic solution. 

is unity, m takes integral values between the limits and we have taken for our 
trial functions M = 2, that is six terms in each component. The half integer 
in the exponent is suggested by the acquisition of a Berry phase after a full 
period -^. For the same reason we have taken the range of the integration 
in equation (2) to be twice the period. So as to create realistic conditions 
for the implementation of the variational procedure, we have chosen a finite 
range for the time variable, although, as can be seen from the formulae in 
equation (26) and equation (27) , for a general value of G the solution is 
not time-periodic. We illustrate the procedure for two cases. 



4.1.1 A periodic case 

This comes about when G is such that the Rabi frequency K in equation 
(27) and u are commensurate. We have chosen G = \/TEuj, so that K = 2uj. 



Minimizing the action with respect to the coefficients subject to the 
normalizing conditions turns out to be equivalent to diagonalizing a twelve- 
by-twelve matrix, whose elements are the action integral computed with the 
circular function shown in equation (28) . There are twelve-by-twelve matrix 
elements, rather than just six-by-six, since the upper and lower states are 
coupled by the off-diagonal terms in equation (25) . 

We find a pair of zero eigenvalues (actually eigenvalues of about 10~^^), 
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Figure 2: Same as in the previous figure but for the imaginary part of the 
amplitudes. 

whose meaning is the value of the action integral in the transformed represen- 
tation; in other words, our variational solution is exact. Also, the numerical 
values of the variationally obtained coefficients A" agree with those in the 
analytic solution in equation (26) . For comparison, the other eigen solutions 
have "eigenvalues" of orders 1-20. 



4.1.2 A non-periodic case 



2^, the analytic solution shown 



With the choice of e.g., G = 3 and K 
in equation (26) is not periodic and an exact solution cannot be achieved 
variationally while having a finite t-range in the action integral. Moreover, 
a larger spread of the basic set in equation (28) is needed. Still, so as to 
estimate the efficacy of the variational procedure under non-optimal condi- 
tions, we used the same t-range and the same set size as in (a). The lowest 
eigenvalues (= the values of the action integral) are about 0.35, compared 
with others eigenvalues, which are again in the range of 1-20. 

The results are also shown graphically, by comparing the variational so- 
lution (full lines) with the exact, analytic solution (broken lines) in equation 
(26) for real and imaginary parts in Figures 1 and 2, respectively. The 
similarities are quite good for the first half (that comprises the period of 
the Hamiltonian) , but is worse in the second half and further deteriorates 
later, say in the time range (— , — ). On the other hand, had we taken the 
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time range of the variational procedure (the upper limit of integration T in 
equation (2) ) up to 77, we would have obviously got a somewhat different 
solution, which would have improved the approximation in the latter range 
and spoiled somewhat the agreement in the earlier range. In general, the 
approximate solution depends on the time-range of the action. In our view, 
this endows a flexibility to the practical application of the method, in the 
sense that, depending on which time range is of interest, better approxima- 
tion can be achieved for that range. Of course, when the "approximate" 
solution is identical with the true solution in some range, it will remain so, 
by analytical continuation, for all ranges. 

In conclusion, one notes a successful application of the variational princi- 
ple for a purely hermitian case, whose solution, though available by algebra, 
is not trivial. 

4.2 A Non-linear Evolution 

The continuous passage of an initially prepared pure state to transitions 
resembling quantum jumps was recently studied in [40], based on a form 
of the Liouville-von Neumann-Lindblad (LvNL) equation. The actual form 
used originated in a representation of fast level crossing in molecular sys- 
tems involving two states [41]. It was noted in [40] that the factorization 
formalism, called the" square root operator" method of [25], represents an 
alternative way of showing how a dissipative term in the Hamiltonian can 
cause decoherence. To apply our variational formalism to this case, we first 
formulate the evolution equations in the factorization scheme and solve the 
resulting equations (this is done (A) below) and, secondly, we obtain an ap- 
proximation to the solution by minimizing the action with respect to some 
parameters appearing in the assumed 7's (this is carried out in (B)). 

4.2.1 Decoherence by the "square root operator" method 

The vector potentials consisting of a Hamiltonian and a dissipative (non- 
Hermitian) part are now written for the two factors (71,72) of the density 
matrix as 

Ai = --Gcos(wt)7i - J72 - ir[7i - |72p(l + ^|7iP)/7i" 

A2 = -Gcos(tot)72- J7i-«r[72-|7iP(l-/u|72p)/72^ (29) 
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One notices the similarity of these expressions with the corresponding for- 
mulae in [41] and [40] (where the interpretation of the terms is spelt out) and 
also the divisor 7 on the extreme right, characteristic of the factorization 
formalism for dissipative processes, equation (10) . The trace of the density 
{p — 77^) stays constant during the motion and this property is maintained 
unchanged also by the super linear terms which enter with the coefficient ^u. 

It may be noted that the equations of motions for 71 (t) and 72 (t) and of 
their conjugates, given in equation (7) , lead to the master equations for the 
matrix elements of the density operator. For the diagonal elements these 
are of the LvNL form (when ^ = 0), but not for the non-diagonal ones. This 
property has already been noticed in [25]. 

We next solve two equations for the 7's, subject to the pure state initial 
conditions 71 = 1, 72 = at t = 0, then form from the solutions the diagonal 
matrix elements of the density matrix and finally show the results in figure 
3. The quantity changed between the upper three drawings is the strength 
r of the dissipative term. As this increases, a transition takes place from 
the slow to the fast decoherence regime. We note the remarkable similarity 
of the results obtained here by the factorization method to those in Figure 
1 of the above articles, except that for strong dissipation their drawings 
show little oscillations, unlike our third drawing from above. In this, drawn 
for — = 20, after a very steep initial slope (not visible in the figure), both 
diagonal density matrix elements oscillate about the asymptotic value of 2- 

The majority of calculations whose results are shown in this section are 
carried out for a density matrix referring to an " assembly" consisting of one 
system. This means that in Appendix A, one has A'^ = 1 and the summation 
over a in equation (34) is trivial. We have also worked the density matrices 
when there is a non-trivial summation, namely when initial conditions on the 

JTT 2i7r SJTT 4i7r 

two factors are 7i(0) =62,62^62^62 and 72(0) = 0, (so that N = 4), 
instead of having only 71(6) = 1 (and A = 1), as before. The resultant 
density tends now to an almost perfectly straight line. This is similar to the 
graphs shown in both [41] and in [40] for strong dissipation and elucidates 
the meaning of system-averaging in the density matrix. (In more complex 
systems, one would require an averaging in the density matrix over a much 
larger number of states, such that A -^ 00.) We have also worked out the 
A = 4 case for the upper three graphs in figure 3. For these graphs, there 
was hardly any perceptible change from those shown. 

For the detailed interpretation of these results, which is not the primary 
subject of this work, we again refer to [40, 41] and references therein. Also, 
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we do not delve here into the details or extensions of the results, but turn 
to the variational treatment of time-development to be got from minimizing 
the action. 

4.2.2 Solving variationally 

We note that the strongly oscillatory factor in the solution arises from the 
driving term ^cos{iot) and that this term was already present in the Hamil- 
tonian case considered in the previous subsection, in equation (26) . (We 
have, however, eliminated there the fast oscillating factor by subtracting 
from the Hamiltonian the so-called dynamic-phase.) So in this subsection 
we shall put G = 0, which also makes the numerical aspect of the variation 
considerably simpler. We then set up pair of suitable trial 7(t)'s, containing 
parameters to be varied. 

In contrast to the previous subsection (which was a linear problem and 
in which a large number of Fourier coefficients were varied), in the present 
problem only one variational parameter (v) is introduced. However, to make 
progress, we must consider critical regions of the time domain, namely t = 
and t = oo. At the former, it is easy to see that in order that the singularity 
due to the zero divisor 7'''2(0) in the vector potential be matched by the 
time derivative of 72 (t) at t = 0, this function must take there the form of 

lim 72(t) -^ Vie'f^Tt + 0(t) (30) 

with the constant phase angle /3 arbitrary. Similarly, it can be shown that, 
asymptotically for large t, the same solution must have the form 



lim 72(t) ^ e-^-^S/[l + t,e-2(r-iJ)i + o(e-2rt)]/2 (31) 

or some other form equivalent to this. A constant phase factor was ignored 
here. The parameter v cannot be found from the equations of motion. We 
seek to obtain the variationally best v, such that the condition at t = is 
also satisfied. After some elementary algebra on finds that in terms of the 
variational parameter, the density factor 72 can have the form. 



72(t) = e-^-^V[l + fe-2(r-*^)t - (l + t,)e-«"+')^-'"-^)— ]/2 (32) 

The previously introduced phase f3 was varied independently and found to 
be small. So we put (3 = 0. The other density factor 71 (t) was so constructed 
that the exponentials inside the square root had the opposite signs to those 
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Figure 3: Unitary evolution. The diagonal element of the density matrix p22 
under a successively increasing dissipative parameter T. The parameters in 
equation (29) are G = 45 (25 in the bottom drawing, chosen to facilitate 
the computation), J = 3, uj = 1, fi (superlinearity coefficient) = 0. Then, 
in top drawing: F = 0, A^ (Ensemble size)= 1; in second drawing: F = 0.05, 
A^ = 1; in third drawing: F = 20, A^ = 1; in bottom drawing F = 20 (as 
in previous, but) A'^ = 4. The initial slope in the last two drawings is too 
steep to be visible and so are tiny fluctuations in the horizontal part of the 
bottom drawing. A non-zero value of the superlinear parameter ^u, of the 
order of 1, hardly changes the curves. 
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in 72 and normalizing factors were added so that |7iP + I72P = 1 at all 
times. 

Minimizing the action integral for a set of parameter values for J and 
r yields optimized f's. Keeping J = 3 fixed (as in references [41, 40] and 
selecting a set of dissipation parameter we obtain as follows: 

J = 3, r= 0.1, V = -.985 

J = 3, r = 0.5, V = - .78 

J = 3, r=2, V = - 0.755 

For the middle case F = .5 we compare in Figure 4. the variational 
solution (thin continuous curve, with v = —.78) and the solution from the 
equation of motion (thick continuous curve). In the asymptotic regime of 
large t, the behavior of the two curves is quite similar, though the amplitude 
of oscillation is clearly smaller in the variational solution than in the exact 
solution. The discrepancy appears the more serious in that a different choice 
of the variational parameter, namely v = —0.98, (also shown in the figure 
by broken lines), which has an action larger than the optimized one, comes 
nearer to the exact solution. However, we show in the inset that in the 
extremely short time region, the optimized solution is qualitatively better 
than the other choice. Because of the singular behavior of the short time 
region in the vector potential, this region dominates the value of the action. 
At the same time, yet another choice of the parameter v = —0.05, (shown 
in the figure by dotted lines) gives a definitely poorer resemblance to the 
correct curve. 

In conclusion, when it comes to describe subtle quantum mechanical en- 
semble properties, the factorization (or "square root operator") method can 
be used either in its equation of motion form or variationally. In the present 
case, at least, the equation of motion approach was from a numerical view- 
point undoubtedly superior. So one may question the use of the variational 
method. However, on the one hand not all problems may be easily solvable. 
On the other hand, one should remember that the derivation of the equation 
of motion in equation (7) is itself based on a variational ansatz introduced 
in this paper. 

5 Further Extensions 

To treat non-Markovian processes, the vector potentials have to be functions 
of the 7 vectors at earlier times, but otherwise, no change in the formalism 
is needed. 
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Probability 




Figure 4: Comparison of solutions. p22) a diagonal component of the 
density matrix, for the following values of the parameters in equation 29: 
G = 0, J = 3, uj = 1, r = 0.5, fi = 0. Solution of equation of motion: 
thick line. Variational solution in equation (32) with optimized parameter 

V = —0.78: thin line. Hypothesized solution as in equation (32) with 

V = —0.98: small broken lines. Hypothesized solution with v = —0.05: large 
broken lines. The inset shows the curves near t = 0. 
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Non-negativity of the entropy-change fohows from the master equations 
and properties of the scattering probabiUties M^a in equation (12) , as is 
shown in [42]. 

Transport processes can be treated simply. Thus, let us consider elec- 
tronic conduction in a solid due to a spatial gradient in the potential (i.e., 
an electric field) or in the ambient temperature. The 7 vectors are, nor- 
mally, labeled by the reciprocal, /c-vector index and are essentially small 
deviations from the square root of the equilibrium, Fermi-Dirac electronic 
distribution function. Following the standard treatment given in, e.g., [13], 
the time derivatives of the 7 vectors (that are now real and identical to 7'^) 
are proportional to the spatial gradients. The vector potentials represent the 
scattering integral. Then either equation of motion in equation (7) is simply 
the Boltzmann equation in an inhomogeneous form; namely, its left hand 
side represents the source or the gradient and the right hand side contains 
the desired distribution function under the integral over all wave- vectors. 
The Lagrangian can be used to obtain the solution variationally. This vari- 
ational formulation is, however, different from those given in [11, 13]. (Of 
course, different variational procedures can lead to the same result.) 

We have noted earlier that the postulated Lagrangian does not contain 
a potential. Adding a potential to the Lagrangian might apparently change 
the equations of motion. It seems to us, however, that under conditions pre- 
vailing in stochastic processes, this will not happen. The reason for this can 
be stated in various forms and is rooted in the circumstance (already noted 
above) that in the presence of a random force one has no control over the 
value of the variables, only on its rate of change. ("Free terminal endpoint" 
condition of [3], p. 9.) In Appendix B we give a formal proof for the follow- 
ing proposition: "When the following conditions hold: (a) the potential is 
a non- negative quadratic form in all of its variables (the 7's), (b) the vector 
potentials A are all real and positive, (c) the initial values of the 7 vectors 
are suitably chosen, and (d) the variation is performed under conditions of 
fixed initial values of 7 and 7, then it follows that the action obtained from 
the variation of the velocities only, i.e., with the potential regarded as ignor- 
able, is less than the action obtained from the variation of both the variables 
and their derivatives (namely, through the usual Lagrange equations, which 
are obtained under fixed initial and final boundary conditions)." The result 
holds probably under a wider range of conditions, since in the proof we have 
not utilized the requirement imposed on the self-correlation of the random 
forces fa by the "second" theorem of Mori [43]. (This requirement ensures, 
among other things, the time-shift invariance of the random process which 
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is at the root of the Onsager-Machlup theory [5, 6]. Needless to say, that 
the result obtained in Appendix A is not in conflict with the validity of the 
Euler-Lagrange equations, since these are obtained under conditions that 
the variables have fixed values at the final time. 

6 Conclusion 

The variational action (or Lagrangian) proposed in equation (2) for dynam- 
ical processes has the advantages of being simple, general and fiexible. It 
differs from previously employed variational procedures by the factorization 
ansatz in equation (1) , by the absence of a scalar potential term and the 
presence of a variable final time upper limit. The relation of the postulated 
Lagrangian to some basic invariance property (like " frame indifference" [44] ) 
remains to be explored, account being taken of the fact that, for vector po- 
tentials that are not all equal, the formalism is non-Abelian (namely, the 
vector potentials cannot be transformed away by a single gauge factor) [45]. 
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A A Tutorial on the Factorized Density Matrix 
Formalism. 

Though the factorized density matrix, written in an abstract form as p = 
7-7"'", has been employed before in [25] and [16], we shall explain its formal- 
ism here, following Band's introductory texts to the von Neumann matrix 
method ([46], [4 7]). Let ^q, be a possible wave function describing the quan- 
tum state of the a'th system in the ensemble (a = 1,2...A^). It can be 
expanded in terms of a set of eigenstates Un as 

^a=Y.^>n (33) 
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As derived in [46] and other texts, the density matrix arises from the en- 
semble average over all systems in the sense that its nm component is 

/'-- = ^E7:7™+ (34) 

a 

The 7"'s are best viewed as row vectors, distinct for each a (or system) 
and the 7"''''s as column vectors. The 7- and 7'''- derivatives in the text 
(which implement the variation procedure) are with respect to 7° and 7^^. 
The ensemble averaging, namely the summing over a and the subsequent 
division by A^, is not explicitly written out in the text, but is designated by 
inserting a dot between 7 symbols, so that the previous matrix element is 
written as 

Pnm = ln-lm (35) 

It is clear that the 7's are not vector quantities, but the traces over the 
dotted products are proper scalars. 

B Proof of the minimal action under one-point 
boundary condition 

Assumptions: In the action, equation (2) , the vector potential A is now as- 
sumed to be positive (non-negative) and real. We shall further subtract from 
the action (see below) a potential term, in which the potential V depends 
on the variables only, not on their derivatives. This potential is supposed 
to be monotonic, non-decreasing and positive in the relevant range of its 
variables. 

We start the proof for a single time dependent variable g which replaces 
the earlier complex variable 7 through 

g = g{t) = -ij (36) 

The reality of g for all times will be evident. The one-point (initial time) 
boundary conditions fix g{0) and ^(0), while g{t) at later times develops 
according to its equation of motion. We write the action, including the 
potential, in the single variable g as 

S{T)= rdt[{g-A{t)f-2V{gm (37) 

Jo 
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The boundary conditions fix the value oi g{0) > and oi g{0). We next min- 
imize the above action in two ways and subsequently compare the resulting 
actions. The first is the usual Lagrange equation way in the presence of a 
potential V and the quantities arising from this method will be denoted by 
the superscript V. The second method pretends that there is no potential 
and the corresponding quantities will take the superscript 0. It is the second 
method that was used in the text. 

ff^(t) = Ait) - V'igit)) (38) 

where the prime represents the derivative with respect to the argument g. 

g\t) = A{t) (39) 

The latter equation imposes the following initial condition for the velocities: 

5^(0) = A(0) 

= ^''(O) (40) 

Integrating equation (38) once, we obtain 



gV^t) = A{t) - P{t) P{t) ^ V'{g{t))dt (41) 



where P{t) is zero at t = and is for positive times non-negative, since it 
is, by equation (38) , the time integral of a positive quantity. Subtracting 
(p shown in equation (39) from the last equation and integrating, it is clear 
that g^ never exceeds (algebraically) g^ . Calculating the actions obtained 
in the two methods and subtracting we find: 

5^(T) = r dt[{g'' - A{t)f - 2V{g^{t))] = f dt[P{tf - 2V{g^m 
Jo Jo 

S\T) = r dt[{g^ - Ait)f - 2Vig\t))] = (^ dt[-2V{g\t))] 
Jo Jo 

=> 5^(T)-50(T)= r dt[Pitf + 2iVig\t)) -Vig^'itm (42) 

Jo 

In the integrand the squared term is necessarily positive (non-negative) and 
so is the term containing the difference of potentials since the 0-argument 
is larger than the V-argument and the potential is monotonic by supposi- 
tion. Though obtained under restricted conditions, the result shows clearly 
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that the two-point boundary conditions are necessary requisites for the va- 
Udity of the Lagrange-Euler equations of motion. Generahzation to several 
(real) variables gi,g2, -.-jgN is immediate, when the potential is a positive 
quadratic form in these variables, since this can be diagonalized (with posi- 
tive eigenvalues) simultaneously with the kinetic energy term. However, the 
initial point variables need to be chosen carefully in this case. 

Finally, we have not proven that the action using equations of motion of 
the text is minimal, but only that is lower than that obtained with the (for 
this case, inappropriate) use of the Lagrange equations. Furthermore, it is 
not evident that the solutions obtained in this Appendix satisfy conditions 
required from density matrices or probabilities (e.g., normalizations). 
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